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SUMMARY 

Seismic studies indicate that the Earth’s inner core has a complex structure and exhibits a strong 
elastic anisotropy with a cylindrical symmetry. Among the various models which have been pro¬ 
posed to explain this anisotropy, one class of models considers the effect of the Lorentz force 
associated with the magnetic field diffused within the inner core. In this paper we extend pre¬ 
vious studies and use analytical calculations and numerical simulations to predict the geometry 
and strength of the flow induced by the poloidal component of the Lorentz force in a neutrally 
or stably stratified growing inner core, exploring also the effect of different types of boundary 
conditions at the inner core boundary (ICB). Unlike previous studies, we show that the bound¬ 
ary condition that is most likely to produce a significant deformation and seismic anisotropy is 
impermeable, with negligible radial flow through the boundary. Exact analytical solutions are 
found in the case of a negligible effect of buoyancy forces in the inner core (neutral stratifica¬ 
tion), while numerical simulations are used to investigate the case of stable stratification. In this 
situation, the flow induced by the Lorentz force is found to be localized in a shear layer below 
the ICB, which thickness depends on the strength of the stratification, but not on the magnetic 
field strength. We obtain scaling laws for the thickness of this layer, as well as for the flow 
velocity and strain rate in this shear layer as a function of the control parameters, which include 
the magnitude of the magnetic field, the strength of the density stratification, the viscosity of 
the inner core, and the growth rate of the inner core. We find that the resulting strain rate is 
probably too small to produce significant texturing unless the inner core viscosity is smaller 
than about 10 12 Pa.s. 
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1 INTRODUCTION 

The existence of structures within the inner core was first discov¬ 
ered by Poupinet et al. (1983), who discussed the possibility of 
lateral heterogeneity from the observation of P-waves travel time 
anomalies. These were then attributed to the existence of seismic 
anisotropy (Morelli et al. 1986; Woodhouse et al. 1986), with P- 
waves travelling faster in the north-south direction than in the equa¬ 
torial plane. Since then, more complexities have been discovered 
in the inner core: a slight tilt in the fast axis of the anisotropy, 
radial variations of the anisotropy with a nearly isotropic up¬ 
per layer, hemispherical variations of the thickness of the upper 
isotropic layer, an innermost inner core with different properties in 
anisotropy or attenuation, and anisotropic attenuation (See Souriau 
et al. 2003; Tkalcic & Kennett 2008; Deguen 2012; Deuss 2014, 
for reviews, and references therein). 

The seismic anisotropy can be explained either by liquid in¬ 
clusions elongated in some specific direction (shape preferred ori¬ 
entation, SPO) (Singh et al. 2000) or by the alignment of the iron 
crystals forming the inner core (lattice preferred orientation. LPO). 
In the case of LPO, the orientation is acquired either during crystal¬ 
lization ( e.g. Karato 1993; Bergman 1997; Brito et al. 2002) or by 
texturing during deformation of the inner core. Several mechanisms 


have been proposed to provide the deformation needed for textur¬ 
ing: solid state convection (Jeanloz & Wenk 1988; Weber & Ma- 
chetel 1992; Buffett 2009; Deguen & Cardin 2011; Cottaar & Buf¬ 
fett 2012; Deguen et al. 2013), or deformation induced by external 
forcing, due to viscous adjustment following preferential growth at 
the equator (Yoshida et al. 1996, 1999; Deguen & Cardin 2009), or 
Lorentz force (Karato 1999; Buffett & Bloxham 2000; Buffett & 
Wenk 2001). 

Thermal convection in the inner core is possible if its cool¬ 
ing rate, related to its growth rate, or radiogenic heating rate is 
large enough to maintain a temperature gradient steeper than the 
isentropic gradient. In other words, the heat loss of the inner core 
must be larger than what would be conducted down the isentrope. 
However, the thermal conductivity of the core has been recently 
reevaluated to values larger than 90W.m~ 1 .K~ 1 at the core man¬ 
tle boundary and in excess of 150W.m _1 .K _1 in the inner core 
(de Koker et al. 2012; Pozzo et al. 2012; Gomi et al. 2013; Pozzo 
et al. 2014), and this makes thermal convection in the inner core un¬ 
likely (Yukutake 1998; Deguen & Cardin 2011; Deguen et al. 2013; 
Labrosse 2014). Inner core translation, that has been proposed to 
explain the hemispherical dichotomy of the inner core (Monnereau 
et al. 2010), results from a convection instability (Alboussiere et al. 
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2010; Deguen et al. 2013; Mizzon & Monnereau 2013) and is there¬ 
fore also difficult to sustain. 

Compositional convection is possible if the partition coeffi¬ 
cient of light elements at the inner core boundary (ICB) decreases 
with time (Deguen & Cardin 2011; Gubbins et al. 2013) or if 
some sort of compositional stratification develops in the outer core 
(Alboussiere et al. 2010; Buffett 2000; Gubbins & Davies 2013; 
Deguen et al. 2013) so that the concentration of the liquid that 
crystallizes decreases with time. However, the combination of both 
thermal and compositional buoyancy does not favor convection in 
the inner core (Labrosse 2014). 

The strong thermal stability of the inner core resulting from its 
high thermal conductivity (Labrosse 2014) is a barrier to any verti¬ 
cal motion and other forcing mechanisms need to work against it. 
This situation has already been considered in the case of deforma¬ 
tion induced by preferential growth in the equatorial belt (Deguen 
& Cardin 2009), and has been shown to produce a layered struc¬ 
ture. Deguen et al. (2011) and Lincot et al. (2014) evaluated the 
predictions of anisotropy from this model and found that although 
it can induce significant deformation, it is difficult to explain the 
strength and geometry of the anisotropy observed in the inner core. 

In this paper, we consider another major external forcing that 
was proposed, Maxwell stress. This was first proposed by Karato 
( 999) who considered the action of the Lorentz force assuming 
the inner core to be neutrally buoyant throughout. This situation is 
rather unlikely and, as discussed above, we expect the inner core 
to be stably stratified. Buffett & Bloxham (2000) have shown that 
in this case the flow is confined in a thin layer at the top of the in¬ 
ner core, similar to the case discussed above for a flow driven by 
preferential growth at the equator. However, the growth of the inner 
core gradually buries the deformed iron and this scenario may still 
produce a texture in the whole inner core. All these previous studies 
considered a fixed inner core size and infinitely fast phase change at 
the ICB. The moving boundary brings an additional advection term 
in the heat balance which can influence the dynamics. In the con¬ 
text of inner core convection Alboussiere et al. (2010) and Deguen 
et al. (2013) have proposed a boundary condition at the ICB that al¬ 
lows for a continuous variation from perfectly permeable boundary 
conditions, that was considered in previous studies, to impermeable 
boundary conditions when the timescale for phase change is large 
compared to that for viscous adjustment of the topography. 

In this paper, we investigate the dynamics of a growing inner 
core subject to electromagnetic forcing, and include the effects of 
a stable stratification, of the growth of the inner core, and different 
types of boundary conditions. We propose a systematic study of the 
dynamics induced by a poloidal Lorentz force in the inner core and 
develop scaling laws to estimate the strain rate of the flow. 

In Section 2 , we develop a set of equations taking into account 
the Lorentz force and a buoyancy force from either thermal or com¬ 
positional origin. Analytical and numerical results are presented in 
Section 3, scaling laws for the maximum velocity and strain rate 
are developed in Sections 4 and 5 and compared to numerical solu¬ 
tions. In Section 6, we use our results to predict the instantaneous 
strain rates and cumulative strain in the Earth’s inner core due to 
the Lorentz force. 

2 GOVERNING EQUATIONS 

2.1 Effect of an imposed external magnetic field 

The magnetic field produced by dynamo action in the liquid outer 
core extends up to the surface of the Earth, but also to the center- 


most part of the core. Considering for example a flow velocity of 
the order of the growth rate of the inner core gives a magnetic 
Reynolds number (comparing advection and diffusion of the mag¬ 
netic field) of the inner core of about 10 -5 . This shows that the 
magnetic field in the inner core is only maintained by diffusion 
from its boundary. 

Two dynamical effects need to be taken into account: the 
Lorentz force and Joule heating. The Lorentz force acts directly 
on the momentum conservation, while Joule heating is part of the 
energy budget and modifies the temperature distribution, inducing 
a flow through buoyancy forces. 

In this paper, we will discuss the effect of the Lorentz force 
in the case of a purely toroidal axisymmetric magnetic field with a 
simple mathematical form. The effect of Joule heating in the case 
of a non growing inner core was studied by Takehiro (2010) and 
will not be investigated further here. 

The poloidal magnetic field intensity at the core mantle bound¬ 
ary (CMB) can be inferred from surface observations of the field 
at the Earth’s surface, but both poloidal and toroidal components 
are poorly known deeper in the core. The root mean square (RMS) 
strength of the field at the ICB has been estimated using numeri¬ 
cal simulations to be around a few milliteslas (e.g. Glatzmaier & 
Roberts 1996; Christensen & Aubert 2006). It can be also con¬ 
strained by physical observations: for example, Root & Dumberry 
(2013) give an upper bound of 9-16mT for the RMS field at the 
ICB by looking at the dissipation in the electromagnetic coupling, 
while Gillet et al. (2010) suggest 2-3mT from the observation of 
fast toroidal oscillations in the core. Buffett (2010) obtains similar 
values from measurements of tidal dissipation. Numerical simula¬ 
tions also predict a strong azimuthal component B ^ at the vicinity 
of the inner core, possibly one order of magnitude higher than the 
vertical component B z (Glatzmaier & Roberts 1996), though this 
depends on the magnitude of inner core differential rotation. 

Buffett & Wenk (2001) have considered the effect of the az¬ 
imuthal component of the Lorentz force resulting from the com¬ 
bination of the B z and B$ components of the magnetic field. We 
will focus here on the effect of the azimuthal component of the 
magnetic field, for which the associated Lorentz force is poloidal 
and axisymmetric. The flow calculated by Buffett & Wenk (2001) 
is decoupled from the flow induced by the azimuthal component 
of the magnetic field, and thus the total axisymmetric flow can be 
obtained by simply summing the two flows. 

One of the most intriguing feature of the Earth’s magnetic field 
is the existence of reversals. However, since the Lorentz force de¬ 
pends quadratically on the magnetic field, its direction is not mod¬ 
ified by a reversal of the field. For simplicity, we will consider that 
the magnetic field is constant in time. 

The magnetic field inside the inner core is calculated by diffus¬ 
ing the field from the ICB. The magnetic Reynolds number for the 
inner core being very small, B is not advected by the flow. Because 
the seismic observation of anisotropy is of large scale, and also be¬ 
cause low-order toroidal component penetrates deeper inside the 
inner core, only the lowest order of the azimuthal component of the 
magnetic field is taken into account, following the work of Karato 
(1999) and Buffett & Bloxham (2000). 

We consider a purely toroidal axisymmetric field of de¬ 
gree two in the vicinity of the ICB, of the form B\ ICB = 
Bo sin # cos (Buffett & Bloxham 2000). Solving V 2 S = 0, 
the field inside the inner core is 

r 2 

B = Bo cos 9 sin Oej,, (1) 

r- 
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Figure 1. Meridional cross sections showing the intensity of the magnetic 
field (a), the Lorentz force F /. (b) and the non-potential part of the Lorentz 
force f i as defined in equation (3) (c). 


in spherical coordinates, which is associated to an electric current 
density J = fr V x B, where r- lc is the radius of the inner core 
and p .o is the magnetic permeability. 

The Lorentz force is a volume force given by Fl = J x B. 
The Lorentz force can be decomposed as the sum of the gradient of 
a magnetic pressure and a non-potential part as Fl = - VT m + 
f L , which is a unique Helmholtz decomposition for V ■ f L = 0. 
With the magnetic field as defined in Eq. (1), we find that P m and 
f L are given by 


IBfrf 
7 Uo rf 



e + 


( 2 ) 


and 


Sl = 


Bo r 3 r 
fiQVic r? 


Q 4 n 7^ 2 n , 4 

3 COS 0 -— COS 0 + — 


+ cos 9 sin 0 ( - — 3 cos 2 6 ) e g 


(3) 


The potential part of the Lorentz force can only promote a 
new equilibrium state but no persisting flow. We are thus only inter¬ 
ested in the non potential part of the Lorentz force, shown in Fig. 1 . 
Eq. (3) provides a characteristic scale for the force as Bq / pori c . 

Karato ( 1 999) investigated the effect of the Maxwell stress by 
applying a given normal stress on the inner core boundary. This is 
different from our study, where, as in Buffett & Bloxham (2000), 
we consider a volumetric forcing, as shown on Fig. 1, and not a 
forcing on the surface of the inner core. 


2.2 Conservation equations 

2.2.1 Conservation of mass, momentum and energy 

We consider an incompressible fluid in a spherical domain, with a 
newtonian rheology of uniform viscosity rf, neglecting inertia. Vol¬ 
ume forces considered here are the buoyancy forces, with density 
variations due to temperature or compositional variations, and the 
Lorentz force as defined above. 

The equations of continuity and conservation of momentum 
are written as 


(4) 


0 = -Vp' + Apg + pV 2 u + f L , (5) 

where u is the velocity, p' the dynamic pressure that also includes 
the magnetic pressure, A p the density difference compared to the 
reference density profile, and g = gi c r/ri c e r the acceleration of 
gravity with <j; c the acceleration of gravity at r — n c . 

The density depends on both the temperature T and the light 
element concentration c. We define a potential temperature as 0 = 
T — T 3 (r,t), with T 3 (r,t) the isentropic temperature profile an¬ 
chored to the liquidus at the ICB, and introduce a potential compo¬ 
sition C = c—c£ c (t), where cf c (t) is the composition of the solid at 
the ICB. We will consider separately the effects of composition and 
temperature, but both can induce a density stratification, which is 
quantified through a variation of density A p which is either par& 
or pacC, where p is the reference density, and ar and ac the 
coefficients of thermal and compositional expansion, respectively. 
Because the potential temperature and composition are solutions of 
mathematically similar equations, we will use a quantity \ which is 
either the potential temperature 0 or composition C. In this paper, 
quantities that apply for both cases will have no subscript, whereas 
we will use t for quantities referring to the thermal stratification, 
and c for compositional stratification. 

The momentum conservation equation (5) is thus written as 

0 = -Vp + apxgtc — e r + tfV 2 u + f L . (6) 

r ic 

The equations for the evolution of the potential temperature (energy 
conservation) and of light element concentration (solute conserva¬ 
tion) have a common form, which will be written as 

|jr + u- Vx = «V 2 x + S(t), (7) 

where k is the diffusivity of either heat (kt) or composition (kc) 
and S a source term given by 

BT 

S T (t) = k t V 2 T s - (8) 

and 

BE 

Sc{t) = (9) 

Bt 

As discussed in Deguen & Cardin (2011), the inner core is 
stably stratified when the source term S(t) is negative, and no con¬ 
vective instability can develop. In this paper, we will focus on this 
case, with either Sr(t) or Sc{t) negative. 


2.2.2 Growth of the inner core 


To take into account the growth of the inner core, we use a front fix¬ 
ing approach to solve the moving boundary problem (Crank 1984) 
by scaling lengths with the inner core radius ri c (t) at time t. We 
define a new coordinate system with r = r/ri c (t). This modifies 
slightly the spatial derivatives by bringing a factor 1 /rt c (t) to radial 
derivatives, but also adds a radial advection term in the equations 
where the time derivative is present. In the new coordinate system, 
we obtain 



B . Uidt ) d 

dt r r r ic (t) dr t 


( 10 ) 


where Ui c (t) = dri c /dt is the instantaneous growth rate of the 
inner core. Eq. (7) becomes 



1 

n c {t) 


(u 


ruic{t)e r ) ■ Vy 


rl(t) 


V 2 X + S(f), 


(ID 


V ■ u = 0, 
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where V and V 2 are now spatial derivative operators in the new 
coordinate system ff, 6 , cf>), with 6 and (j> the colatitude and longi¬ 
tude. 


2.3 Dimensionless equations and parameters 

2.3.1 Definition of the dimensionless quantities 

The set of equations (4), ( 6 ), (7) is now made dimensionless, us- 
ing Vic ft), the age of the inner core n c , n/ricft), V K / r icft) and 
A p x as characteristic scales for length, time, velocities, pressure, 
and density variations. The density scale A p x is the difference of 
density across the inner core due to either thermal or compositional 
stratification. The quantity x is scaled by A p x /ap. The character¬ 
istic velocity scale is defined using the diffusion time scale rather 
than the inner core growth rate, to make it usable in both the grow¬ 
ing and non-growing inner core cases. The quantity S(t) is made 
dimensionless using rf c ap/nAp x . Using the same symbols for the 
dimensionless quantities (including using now r for the dimension¬ 
less radius f defined in the last subsection), we obtain 

V ■ u = 0, (12) 


0 = -Vp' + Ra(t)xre r + V 2 u + Mft)f L , (13) 


m 


dx 

dt 


(u — Pe{t .) r e r ) • Vx + V 2 x 


+ Bft) 



(14) 


with four dimensionless parameters defined as 


Raft) = 

Apx (t)gicfic(t) 

T]K, 

(15) 

Mft) = 

Bo r ic{t) 

(16) 

POVK 

m = 

r? c (t) 

KTic 

(17) 

Pe(t) = 

UicT'ic(t ) 

K 

(18) 


The last term in Eq. (14) comes from the time evolution of the scale 
for x, A p x /ap. 

£(t) and Pe(t) characterize the growth of the inner core. The 
Peclet number Pe(t) compares the apparent advection from the 
moving boundary to diffusion. A large Peclet number thus corre¬ 
sponds to a fast inner core growth compared to the diffusion rate. 
In the case of a non-growing inner core, Pe = 0, S(t) = 0 and 
the relevant timescale is no longer Ti c but the diffusion time scale, 
which gives £ = 1. This approach allows us to treat both non¬ 
growing and growing cases with the same set of dimensionless pa¬ 
rameters. 

M ft) is an effective Hartmann number, which quantifies the 
ratio of the Lorentz force to the viscous force, using thermal dif- 
fusivity in the velocity scale. This effective Hartmann number is 
related to the Hartmann number often used in magnetohydrody¬ 
namics (Roberts 2007), Ha = Br/goqX, through M = Ha 2 A/«, 
where A is the magnetic diffusivity. 

Raff) defined in equation (15) is the Rayleigh number that 
characterizes the stratification, and is negative since A p x is neg¬ 
ative for a stable stratification. The density stratification depends 
on the importance of diffusion, and on the time-dependence of the 


inner core radius. Expressions for Apr and A p c will be given in 
section 2.3.2. 

To solve numerically the momentum equation (13), the ve¬ 
locity field is decomposed into poloidal and toroidal components. 
The complete treatment of this equation and the expression of the 
Lorentz force in term of poloidal and toroidal decomposition are 
described in appendix A. 


2.3.2 Simplified growth of the inner core 

A realistic model for the inner core thermal evolution can be ob¬ 
tained by calculating the time evolution of the source term St ft) 
and the radius ri c ft) from the core energy balance (Labrosse 2003, 
2015), as done by Deguen & Cardin (201 1). The result is sensi¬ 
tive to the physical properties of the core. To focus on the effect 
of the Lorentz forces, we choose a simpler growth scenario and 
assume that the inner core radius increases as the square root of 
time (Buffett et al. 1992). Using ry c (f) = r; c (Ti c )(fwith 
fidjic) the present radius of the inner core, the growth rate is thus 

Uicft) — tic. (ti,-) l'2\/ I ■ 

This leads to the following expressions for the control param- 


eters: 



Raff) 

= Ra 0 ^ Px{t) t 3/2 , 

Ap x ,o 

(19) 

Mft) 

= M 0 t, 

( 20 ) 

m 

= 2 Pe 0 1, 

( 21 ) 

Peft) 

= Pe 0 , 

( 22 ) 


where the subscript 0 corresponds to values for the present inner 
core, and t is dimensionless. 

The Peclet number Pe(t) is constant and equal to Peo = 
ficCic)/(2KTic), and the parameter £ is proportional to Pe o. We 
are left with only three independent dimensionless parameters: the 
Rayleigh number Rao characterizes the density stratification, the 
effective Hartman number Mo the strength of the magnetic field, 
and the Peclet number Pe o the the relative importance of advec¬ 
tion from the growth of the inner core and diffusion. 

The value and time dependence of A p x (t) depends on 
whether a stratification of thermal or compositional origin is con¬ 
sidered: 


• In the thermal case, the source term for thermal stratification 
St ft) defined in Eq. ( 8 ) can also be written as 


S T (t) 


pg'^T 

K s 



n c (t)uic(i) - 3 k t , 


(23) 


where dT s / dT a d is the ratio of the Clapeyron slope to the adiabat 
gradient, g' = dg/dr = <?; c /r; c , 7 the Gruneisen parameter, and 
Ks the isentropic bulk modulus (Deguen & Cardin 2011). With 
ri c oc i 1//2 , the product ri c (t)ui c ft) is constant, and so is St- 
Solving the energy conservation equation for the potential temper¬ 
ature (x = 0) assuming u — 0, ri c oc t * 1 ^ 2 , and St constant gives 


0 = 


Srr 2 c 


6kt(1 + Psto/3) 


1 - 


Tie ft) 


(24) 


in dimensional form (see appendix B for the derivation). If Pe 0 <g 

1, then the potential temperature difference A0 across the inner 

core is Srr 2 c /6K, which corresponds to a balance between ef¬ 
fective heating (St) and diffusion. In contrast, A0 tends toward 
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SrTic if diffusion is ineffective and Pe o 3> 1. From Eq. (24), we 
obtain 


Ap T = _ a TpSrrl _ 

6kt( 1 + -PeTo/3) 

and 


(25) 


n _ UTpgi C S T r h ic 

aT + Pe T o/3)' 

With gi c oc ri c and ri c oc t 1 ^ 2 , this ; 


(26) 


.A , 


Rclt ( f ) = RcitoE . 


(27) 


• We estimate the density stratification due to composition from 
the equation of solute conservation, assuming that the outer core is 
well-mixed and that the partition coefficient is constant. The com¬ 
positional Peclet number is large ( Pec ~ 10° with a diffusivity 
Kc ~ 10 -1 ° m.s -2 ) and solute diffusion in the inner core is there¬ 
fore neglected. 

The composition of the solid crystallized at time t at the ICB is 
estimated as 


cLb(f) = kc l 0 



rift) 

r c 



(28) 


(see Appendix C), from which the density difference across the in¬ 
ner core is 


A pc it) = oi C p [cf cb (t) - c 3 icb (t = 0)]. 


= acpkc 0 


1 - 


fic(f) 


- 1 


(29) 

(30) 


We take advantage of the smallness of (ri c (f)/r c ) 3 < 4.3 10 2 to 
approximate A pc as 

3 

(31) 


Apc(t) ~ acpk(l - k)c Q ( Tlc ^ 


which gives 

A p c{t)g ic r 3 c (t) = acpk( 1 - k)c l 0 g ic rf c (t) 
r/Kc 


Rac = 


rjKcrl 


(32) 


With gi c oc ric and ri c oc t 1 ^ 2 , this gives Rac oc rj c oc t 7 ^ 2 and 
Rac it) = Racot 7 ^ 2 . (33) 


2.4 Boundary conditions 


The Earth’s inner core boundary is defined by the coexistence 
of solid and liquid iron, at the temperature of the liquidus for 
the given pressure and composition. By construction, the poten¬ 
tial temperature 0 and concentration C are both 0 at the ICB : 
&{ric(t)) = C(ri c {t)) = 0. The mechanical boundary conditions 
are tangential stress-free conditions and continuity of the normal 
stress at the inner core boundary. 

When allowing for phase change at the ICB, the condition of 
continuity of the normal stress gives 
du 

—V(t)(ur ~ Uic) - 2 -jA- + p = 0 (34) 

in dimensionless form. The parameter P{t) was introduced by 
Deguen et a (2013) to characterize the phase change, and is the 
ratio of the phase change timescale to the viscous relaxation 
timescale t v = r//(Sp g ic bri c ). 


Pit) 


T(j,8p gicbric 


(35) 


where <5p is the density difference between liquid and solid iron at 
the inner core boundary, has been estimated to be ~ 10 3 years 
(Alboussiere et al. 2010; Deguen et al. 2013). The limit V —> 0 
corresponds to perfectly permeable boundary conditions where the 
phase change occurs instantaneously, and V —> oo corresponds to 
perfectly impermeable boundary conditions with no phase change 
allowed at the boundary. 

With ric{t) oc t 1 / 2 and constant, the parameter V{t) is 
expressed using the current value Vo =V[t = n c ) as 

V{t) = Pot. (36) 

2.5 Numerical modelling 

The code is an extension of the one used in Deguen et a (2013), 
adding the effect of the magnetic forcing as in Eq. ( 6 ). The sys¬ 
tem of equations derived in appendix A in term of poloidal/toroidal 
decomposition is solved in axisymmetric geometry, using a spher¬ 
ical harmonic expansion for the horizontal dependence and a fi¬ 
nite difference scheme in the radial direction. The nonlinear part of 
the advection term in the temperature (or composition) equation is 
evaluated in the physical space at each time step. A semi-implicit 
Crank-Nicholson scheme is implemented for the time evolution of 
the linear terms and an Adams-Bashforth procedure is used for the 
nonlinear advection term in the heat equation. 

The boundary conditions are the same as in Deguen et al. 
(2013), but for most of the runs we use P = 10 6 , which correspond 
to impermeable boundary conditions as discussed in section 2.4. 

When keeping the inner core radius constant, the code is run 
until a steady state is reached. Otherwise, the code is run from t = 
0.01 tof = 1 . 

3 FLOW DESCRIPTION 
3.1 Neutral stratification 

In this subsection, we investigate the effect of the boundary condi¬ 
tions on the geometry and strength of the flow by solving analyt¬ 
ically the set of equations in the case of neutral stratification. The 
analytical solution for neutral stratification has also been used to 
benchmark the code for Ra = 0. 

In the case of neutral stratification, with Ra = 0, the equations 
for the temperature or composition perturbation (14) and momen¬ 
tum conservation (37) are no longer coupled. The diffusivity is no 
longer relevant and the problem does not depend on the Peclet num¬ 
ber. Eq. (37) is solved in appendix D using the boundary conditions 
presented in the previous section. The flow velocity is found to be 
proportional to the effective Hartman number M times a sigmoid 
function of P. The solution is shown on Fig. 2, with dimension¬ 
less maximum horizontal velocity and root mean square velocity 
as functions of the phase change number P, as well as streamlines 
corresponding to the two extreme cases, P = 0 (fully permeable 
boundary conditions) and P oo (impermeable boundary condi¬ 
tions). 

In the limit P 1, corresponding to permeable boundary 
conditions, the streamlines of the flow cross the ICB, which indi¬ 
cates significant melting and freezing at the ICB. In contrast, the 
streamlines in the limit P>1 are closed lines which do not cross 
the ICB, which indicates negligible melting or freezing at the ICB. 
The velocity is proportional to the effective Hartmann number M 
whereas the P dependence is more complex. The velocities reaches 
two asymptotic values for low and large P values, separated by a 


V 
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Table 1. Typical values for the parameters used in the text, and typical range of values when useful. 


Parameter 

Symbol 

Typical value 

Typical range 

Magnetic field 

Bo 

3 x 10~ 3 T 

o 

l 

1 

o 

1 

CO 

H 

Thermal diffusivity c 

K'j’ 

1.7 x 10 —5 m 2 .s 

0.33 - 2.7 x 10“ 5 m 2 .s 

Chemical diffusivity a 

K X 

10 -1 ° m 2 .s 

O 

1 

O 

1 

o 

1 

to 

3 

to 

in 

Viscosity 

V 

10 16 Pa.s 

10 12 - 10 21 Pa.s 

Age of IC 

'Tic 

0.5 Gyrs 

0.2- 1.5 Gyrs 

Density stratification (thermal case) d 

Apr 

6kg.m~ 3 

0.5 - 25 kg.m~ 3 

Density stratification (compositional case) 6 

Ape 

5kg.m -3 

1 - 10kg.m -3 

Phase change timescale 

T 0 

10 3 yrs 

10 2 -10 4 yrs 

Inner core radius b 

ric(n c ) 

1221 km 


Acceleration of gravity (r = ri c 

9ic 

4.4m.s -2 


Density of the solid phase b 

p 

12800 kg.m -3 


Density difference at the ICB 

&Pic 

600kg.m~ 3 


Thermal expansivity 

a 

o 

1 

Crt 

1 


Permeability 

Po 

4t r X 10~ 7 H.m -1 



° from Gubbins et al. (2013) 
b from PREM Dzieworiski & Anderson (1981) 

c obtained using k = 163 W. m 1 . K 1 , c p = 750 J.K —1 .kg —1 (Pozzo et al. 2012; Gonti et al. 2013) 
d assuming S = 10 — 1000K.Gyrs _1 (Deguen & Cardin 2011) 
e from Deguen & Cardin (201 1) 


sharp kink. The discontinuity in the derivative of the maximum hor¬ 
izontal velocity slightly above V ~ 10 2 corresponds to a change of 
the spatial position of the maximum, when the streamlines become 
closed and the maximal horizontal velocity is obtained at the top 
of the cell and no longer at its bottom. The change of behavior of 
the boundary from permeable to impermeable induces a significant 
decrease of the strength of the flow, since the velocity magnitude 
in the P>1 regime is one order of magnitude smaller than when 
permeable boundary conditions (V -C 1) are assumed. 

Fig. 2.b shows the maximum of the velocity, now given in 
m.s -1 , as a function of the viscosity, using typical values of the 
parameters given in Table 1 and five different values for the phase 
change timescale r^, from zero to infinite. In term of dimensionless 
parameters, a high viscosity corresponds to small Hartmann num¬ 
ber M and phase change number V. Changing the timescale 
translates the position of the transition between the two regimes, the 
viscosity value corresponding to the transition being proportional 
to Ttj,, but does not change the general trend of the curve, which is 
a linear decrease of the velocity magnitude in log-log space, except 
for the kink between the two regimes. The linear decrease is due to 
the viscosity dependence of the Hartmann number M oc p -1 . For 
typical values of the phase change timescale between 100 years and 
10 000 years, the kink between the two regimes occurs at a viscos¬ 
ity in the range 10 17 — 10 21 Pa.s. 

In what follows, we will focus on the conditions which are the 
most favorable to deformation due to the poloidal component of 
the Lorentz force, and therefore focus on the case of low viscosity 
and large V. The ICB would act as a permeable boundary only if 
V < 10 2 (see Fig. 2), corresponding to rj > 10 1 ' Pa.s. Under these 
conditions, the typical flow velocity would be < 10 _12 m.s _1 , i.e. 
two orders of magnitude or more smaller than the inner core growth 
rate, and would be unlikely to result in significant texturing. For this 
reason, we will let aside the high viscosity/low V regimes to focus 
on low viscosity/high V cases, for which the ICB acts as an im¬ 
permeable boundary. This gives boundary conditions very different 
from previous studies, where perfectly permeable boundary con¬ 
ditions were assumed (Karato 1999; Buffett & Bloxham 2000). In 


particular, this implies that the flow velocity estimated by (Karato 
1999) was overestimated by one order of magnitude. 

According to Eq. (36), the parameter V varies linearly with 
time, which means that V must have been small early in inner 
core’s history. However, this is true for a very short time, when 
the inner core radius was very small, of the order r; c (T; c )/'Pd , 
and this episode is unlikely to have observable consequences in the 
present structure of the inner core. 

3.2 Zero growth rate 

We first investigate the effect of the Lorentz force without taking 
into account the secular growth of the inner core (Pe = 0). Fig. 3 
shows the vorticity and temperature fields obtained for different 
values of the Rayleigh number, at a given effective Hartmann num¬ 
ber M = 10 4 , for a thermally stratified inner core. 

When the Rayleigh number is small, the vorticity field is or¬ 
ganized in two symmetric tores wrapped around the N-S axis.The 
stratification is too weak to alter the flow induced by the Lorentz 
force and the temperature field is advected and mixed by the flow. 
The velocity field is equal to the one calculated analytically for 
Ra = 0 (see subsection 3.1 and appendix D). 

However, when the Rayleigh number is larger, the flow is al¬ 
tered by the stratification and is confined in an uppermost layer, as 
found by Buffett & Bloxham (2000). The velocity is smaller than in 
the case of neutral stratification. The temperature field is strongly 
stratified and the perturbations due to radial advection are small. 
The flow obtained here is similar to the flow induced by differential 
inner core growth with a stable stratification (Deguen et al. 2011), 
with a notable difference: we impose a large V implying a near 
zero radial flow v r across the ICB, whereas Deguen et al. (20 1 ) 
impose a given v r as the driving force. The confinement of the flow 
in a thin layer is likely to concentrate the deformation and thus we 
may expect higher strain rates for a highly stratified inner core, but 
a different spatial distribution of the deformation. 

To explore the parameter space in terms of Rayleigh and effec¬ 
tive Hartmann numbers, we computed runs with Rayleigh numbers 
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Table 2. Typical values of the dimensionless parameters discussed in the text for the present inner core, using typical values from table 1 . 


Dimensionless parameter Symbol Thermal 


Compositional 


Rayleigh number Ra 

effective Hartmann number M 

Peclet number Pe 

Phase change number V 


10 16 Pa.s 


X (-2.8 X 10 s ) 
B 0 \ 2 10 18 Pa.s 


3 X 10- 3 T 

2.8 

10 16 Pa.s 


X 10 4 


10 18 Pa.s 


V 


x 63 


x (-8 x 10 12 ) 

B 0 \ 2 10 18 Pa.s 


3 x 10-3 T 
4.7 x 10 5 
10 18 Pa.s 


x 1.07 x 10 7 


x 10 4 


V 


See the definitions of the dimensionless parameters in the text. 


a. 



Figure 2. Analytic solution for Ra = 0. (a) Evolution of the dimensionless 
velocity as a function of the phase change number V, with streamlines for 
V —> 0 (left) and V —> oo (right). The RMS velocity and the maximum 
of the horizontal velocity are plotted, (b) Evolution of the RMS velocity as 
a function of ij, with velocity in m.s" 1 . Except for the viscosity and the 
phase change time scale the parameters used for definition of V and M 
are given in Table 1. The kink in the curves corresponds to the change in 
regime between large V (low viscosity) and low V (large viscosity), and the 
corresponding viscosity value is a function of the phase change timescale 
T <i>- 


from —10 3 to —10 7 and effective Hartman number from 10° to 
10 6 . The maximum velocity (normalized by M) is plotted in Fig. 4 
as a proxy to determine the regime. The largest velocity coincides 
with the flow velocity obtained for neutral stratification. The vortic- 
ity field corresponding to some of the points in the regime diagram 
are also shown in Fig. 4. 

The systematic exploration of the parameter space reveals 
two different dynamical regimes, which domains of existence in 


Ra = -10 6 Ra = -10 4 



Ra = —10 5 Ra =-6.10 4 Ra = -10 5 Ra = - 6.10 4 
AT = 0.955 AT = 0.277 


Temperature 


Vorticity 

M = 10 4 , Pe = 0 


Figure 3. Snapshots of meridional cross section of the temperature and the 
vorticity fields for M = 10 4 and a constant inner core radius, for four 
different values of the Rayleigh number (from top right, going clockwise: 
Ra = — 10 4 , —6 x 10 4 , —10 5 , —10 6 ). When the stratification is large 
enough (Ra = — 10 6 ), the flow is confined at the top of the inner core and 
the temperature field has a spherical symmetry. When the stratification is 
weak (Ra = —10 4 ), the flow is similar to the one in Fig. 2 for Ra = 0 
and the temperature is almost uniform. The vorticity is scaled by k j- /r? 
and the temperature by Sr? c /6K,T- For u ic = 0, Sv? c /6k,t reduces to 
T s ( 0) - T s (r ic ). 



Figure 4. Maximum velocity (normalized by the Hartmann number M ) in 
the upper vorticity layer, for a zero growth rate. The velocity is scaled by the 
diffusion velocity ft/ri c . The maximum size of the dots corresponds to the 
value for Ra = 0 computed analytically. For some values of (M, —Ra), 
the vorticity field is plotted in the meridional cross section. The red line 
with a slope of 1 shows the limit between the two regimes. 
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a (— Ra,M ) space are shown in Fig. 4. In the upper left part of 
the diagram (large effective Hartmann number, low Rayleigh num¬ 
ber), the flow is very similar (qualitatively and quantitatively) to 
the analytical solution for a neutral stratification, and deformation 
extends deep in the inner core. This regime is characterized by a 
negligible effect of the buoyancy forces, and will therefore be re¬ 
ferred to as the weakly stratified regime. In the lower right part, 
the flow is confined in a shallow layer which thickness depends on 
the Rayleigh number only (not on M ) and in which the velocity is 
smaller than for neutral stratification. This regime will be referred 
to as the strongly stratified regime. 


3.3 Growing inner core 

To investigate the effect of inner core growth, we compute several 
runs with a given set of parameters (Rao, Mo, Pe o), with the time 
t between t = 0.01 and t = 1. Unlike in subsections 3.1 and 
3.2, the dimensionless numbers evolve with time, as described by 
Eqs (19) to (22). 

Fig. 5 shows the evolution of the vorticity field in six sim¬ 
ulations, for a thermal stratification, with the same Rayleigh and 
effective Hartman numbers, Rar o = —10 6 , Mto = 10 4 , but dif¬ 
ferent values of the Peclet number, which corresponds to increasing 
diffusivity from left to right. For each run, snapshots of the vortic¬ 
ity field corresponding to four time steps are shown, from top-right 
and going clockwise. 

Fig. 5 shows that the thickness of the upper layer increases 
with increasing Peclet number. The transition between the two 
regimes of strong and weak stratification is shifted toward larger 
Rayleigh numbers when the Peclet number is increased. At low or 
moderate Peclet numbers (Peo < 10 2 in the cases presented here), 
the magnitude of vorticity is almost constant time, implying that 
the deformation rate in the uppermost layer is also constant. 


4 SCALING LAWS 

In this section, we determine scaling laws for the thickness of the 
shallow shear layer and the maximum velocity in the layer in the 
strongly stratified regime from the set of equations developed in 
section 2. We will first discuss the transition between the strongly 
stratified and weakly stratified regimes discussed in Fig. 4. We will 
then focus on the strongly stratified regime and estimate the defor¬ 
mation in the uppermost layer. Thermal and compositional strati¬ 
fication are discussed separately. The flow in the weakly stratified 
regime is given by the analytical model discussed in section 3.1 and 
appendix D for neutral stratification. 

4.1 Balance between magnetic forcing and stratification 

We start here by discussing the transition between the strongly- 
stratified and weakly-stratified regimes. We base our analysis on 
the vorticity equation obtained by taking the curl of the momentum 
conservation equation (Eq. (13)), 

0 = -Raft f^ge^ + M(t)V x f L + V 2 u>, (37) 

where cj = V x u is the vorticity. The quantity y (denoting ei¬ 
ther potential temperature or composition) is split into two parts, 
y = y(r,t ) + y'(r,9,t), where y is the reference radial profile 
corresponding to u = 0. The vorticity equation then writes 

0 = Ra X e q) + M V x f L + V ut. (38) 

o v 


In the vorticity equation, the three terms must balance if the effect 
of stratification is important. Starting from a state with no pertur¬ 
bations, x' — 0, the flow velocity is initially set by a balance be¬ 
tween the Lorentz force and the viscosity force. Isosurfaces of x are 
deformed by the resulting flow, and the buoyancy force increases, 
eventually balancing the magnetic force if the stratification is strong 
enough. In this case, further radial motion is prevented and the flow 
tends to be localized in a layer below the ICB, as found in our nu¬ 
merical simulations. Denoting by S the thickness of the shear layer 
and u and w the horizontal and vertical velocity, respectively, the 
vorticity is uj ~ u/5. We therefore have, from Eq. (38), 

(-Ra)y'~M~^. (39) 


The perturbation y' thus scales as 

, M 
X ~ —Ra 


(40) 


if the stratification is strong enough for the induced buoyancy 
forces to balance the Lorentz force. 

The effect of the stratification is negligible if the buoyancy 
forces, which are ~ —Ray', cannot balance the Lorentz force, 
which is ~ M. Since y' is necessarily smaller than |x( r ic) — x(0)k 
which by construction is equal to 1, the effect of the stratification 
will therefore be negligible if M —Ra. This is consistent with 
the boundary between the two regimes found from our numerical 
calculations, as shown in Fig. 4, as well as with the results of Buf¬ 
fett & Bloxham (2000) who found that the Lorentz force can dis¬ 
place isodensity surfaces by ~ r- IC M / (—Ra). This estimate is valid 
for both a growing or non-growing inner core. 


4.2 Scaling laws in the strongly stratified regime 

The reference profile y is solution of 


= V 2 y + Per -Vx + S(t)-((^-y. (41) 

Subtracting Eqs (41) to (14), and assuming that y 1 -C y, we obtain 



~fx' 



~ x ‘/s 2 


dy’ dy' 
U d6 W dr 


ux'/S' 





+ Pe r ■ V y — £ 
'-*-' 



~Pex 7 /5 


~Pex' 


(42) 


Three of these terms depend on the growth rate: £dy'/dt, 
Pe r ■ Vy 1 , and £ Ap/Apy'. With our assumption of n c ~ f 1 , 
we have £ = 2 Pe t and thus £ < Pe. Thus, the largest term among 
the growth rate-dependent terms is Pe r ■ Vy' ~ Pe y' /S. 

Comparing the effect of the diffusion term, which is ~ x'/S 2 , 
with the inner core growth term, which is ~ Pey /5, we find that 
the effect of the inner core growth is negligible if 

(43) 

o 

This suggest the existence of two different regimes depending on 
whether Pe is small or large. We develop below scaling laws for 
these two cases. 
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Figure 5. Snapshots of the vorticity field for simulations with dimensionless parameters Mo = 10 4 , Rao = —10 6 , and Pe o = 0.01, 1,10,10 2 ,10 3 and 10 4 
(from left to right), with ri c oc t 1 / 2 . Each panel corresponds to one simulation, with four time steps represented : t = 0.25, 0.50, 0.75 and 1 dimensionless 
time, from top-right and going clockwise. See Fig. 8 for strain rates of corresponding runs. 


4.2.1 Small Pe limit 


Neglecting the growth terms, we have 




06 


''x'7'5 2 


dr 
wx'P 


dr 


(44) 


The conservation of mass implies that u ~ w/8, and with x' ~ 
M/Ra, we obtain 


0 = 



/ RaS 2 



Ox dx 

- w— - w-=- . 

Or or 


r^uM/Ra uM / Ra <~^u5 


(45) 


We now assume that the advection of the perturbation x! is small 
compared to the vertical advection of the reference state, which re¬ 
quires that 5 M/{—Ra). Balancing the advection and diffusion 

terms, we obtain 


M 


u8. 


(- Ra)8 2 

Combining this expression with the relation u/S 3 
from the vorticity equation (Eq. (39)), we have 

c5~ {-Ray 1 / 6 , 

u~ M{-Ra)~ 1/2 . 


(46) 
M obtained 

(47) 

(48) 


4.2.2 Large Pe limit 


In this limit, the diffusion time is larger than the age of the inner 
core, which allows us to neglect the diffusion term. Keeping only 
the largest growth rate-dependent term, and using Eq. (39), we have 


0 = -'‘W-"’W- w W + £ll^L- 

PeM/SRa 


(49) 


Assuming again that the advection of the perturbation is small com¬ 
pared to the vertical advection of the reference state, the main bal¬ 
ance is between the second and third terms, which gives 


u 5 


PeM 
8 {-Ra )' 


(50) 


Combined with u8 ~ MS 4 from Eq. (39), we find that the thick¬ 
ness and maximum velocity of the upper layer are 


M 


Pe 
—Ra ) 

Pe 


1/5 


3/5 


(51) 

(52) 


Two conditions have to be fulfilled for these scaling laws to be 
valid. First, we must have Pe <C 1/8 , which is Pe {—Ra) 1 / 6 
using Eq. (51). Also, we have assumed that the upper layer is thin 
(8 <C 1) and that the horizontal advection is small compared to the 
vertical one. 


4.2.3 Time-dependence 

Our derivation does not make any assumption on the form of the 
inner core growth ri c {t), and the scaling law validity should not 
be restricted to the ty c (i) oc t 1 / 2 case assumed in the numerical 
simulations. 

These scaling laws are valid at all time during the growth of 
the inner core, provided that Pe < 1 /8 and { < 1 and that the 
time dependence of the control parameters is properly taken into 
account. For example, under the assumption of r; c oc t 1 / 2 and with 
M, Ra and Pe given by Eqs. (20), (27) and (22), we obtain 

8 ~ {-Ra T0 r 1/6 1~ 1/2 , u ~ MTo{-RaTo)~ 1,2 t~ 1/2 (53) 

for the thermal case (small Per), and 


8 


f Peco \ ^ ,- 7 /io 

\{-Ra C o)J 


u 


Mco 


Peco 

{—Raco) 


t -n/io 

(54) 


for the compositional case (large Pec). 


4.3 Comparison with numerical results 

Fig. 6 shows the thickness of the uppermost vorticity layer and 
the maximum horizontal velocity obtained in numerical simula¬ 
tions with a constant inner core radius, which corresponds to the 
small Peclet number limit. When —Ra/M -C 1, the flow has 
the geometry and amplitude predicted by our analytical model for 
Ra = 0. When —Ra/M > 10 and the thickness of the upper layer 
is smaller than ~ 0.3, the data points align on straight lines in log- 
log scale, with slopes close to the predictions of the scaling laws 
developed in subsection 4.2.1 (Eqs. (47) and (48)). 

Fig. 7 shows the vorticity layer thickness and maximum ve¬ 
locity as functions of Pe/{—Ra), in log-log scale, for runs in the 
large Peclet limit. The thickness of the upper layer and the maxi¬ 
mum velocity align on slopes close to the 1/5 and 3/5 slopes pre¬ 
dicted in subsection 4.2.2. Fig. 7 has been constructed from runs 
with M = 1, but we have checked that, as long as the condition 
—Ra M is verified, the geometry of the flow does not depend 
on M and that the velocity is proportional to M. 


u 


—Ra 
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Figure 6. Results from simulations with a constant inner core radius. Evolu¬ 
tion of the thickness of the uppermost vorticity layer (a) and maximum hor¬ 
izontal velocity (b) with the absolute value of the Rayleigh number —Ra. 
Colors correspond to different values of the effective Hartmann number M, 
and dashed lines to the corresponding —Ra = 10 M line. The velocity has 
been scaled with the effective Hartman number and the extreme value for 
low —Ra corresponds to the analytical model with no stratification (black 
horizontal dotted line). The solid black lines are the best fit for 5 < 0.3, 
5 = 1.81(-Ra)-° 16±0 013 and V max = 0.44M(-_Ra)°- 506± ° 002 . 
The orange lines show the slopes predicted in subsection 4.2.1. 


a. 



Figure 7. Thickness of the uppermost vorticity layer (a) and maximum ve¬ 
locity (b) as functions of —Pe/Ra. 50 runs are plotted, with Rao from 
—3 x 10 3 to —10 10 and Peo from 13 to 5000, with 10 time steps for 
each runs, and filtered by Peo 1 /S. Solid green lines are the best fit 
for 5 < 0.25, which is <5 = 1.24 (-Pe/Ra )°■ 192±0 - 04 and V mdX = 
0.143 M (— Pe/Ra)°- 56±0 09 . The red lines are expected scaling laws 
developed in subsection 4.2.2. 


5 STRAIN RATE 


Fig. 8 shows the von Mises equivalent strain rate for the runs corre¬ 
sponding to Fig. 5, highlighting regions of high deformation. The 
von Mises equivalent strain rate is the second invariant of the strain 
rate tensor, measuring the power dissipated by deformation (Tome 
et al. 1984; Wenk et al. 2000; Tackley 2000). Comparing Figs 5 
and 8, we see that the deformation and vorticity fields have a simi¬ 
lar geometry when the flow is organized in several layers, whereas 
the location of the regions of high deformation and high vorticity 
differ when the effect of stratification is small and the flow is orga¬ 
nized in one cell only. In this case, which is similar to that studied 
by Karato (1 999), the maximum deformation is at the edges of the 
cells, whereas in the case of large stratification, the strain is con¬ 
fined in the uppermost layer. In the strongly stratified regime, the 
deformation can be predicted from the scaling laws discussed in 
section 4.2, as i ~ u/8 in dimensionless form, or i ~ un/8r 2 c in 
dimensional form. 

In the small Pe number case, relevant for the Earth’s inner 
core with a thermal stratification, using the scaling laws (47) and 
(48) for the velocity and shear layer thickness gives 


e(t) 


K, U 


~ 0.2^M(f)(-i?a(f))- 1/3 . 


(55) 


In the large Pe number case, relevant for the Earth’s inner core 
with a compositional stratification, using the scaling laws (5 1 ) and 
(52) for the velocity and shear layer thickness gives 


e(t) 


K U 




/ —Ra(t) 

V Pe(t) 


-2/5 


(56) 


Notice that i is proportional to r /~ 2 ^ 3 and uj ~ 3 ^ 3 in the thermal 
and compositional cases, and therefore increases with decreasing 
viscosity in spite of the fact that the velocity in the shear layer de¬ 
creases with decreasing r/. This is because the thickness of the shear 
layer decreases with viscosity faster than the velocity. 

These scaling laws give upper bounds for the actual strain rate 
in the inner core, which evolves with time because of the time de¬ 
pendence of the parameters. The quantity 1/e is the time needed to 
deform the shear layer to a cumulated strain ~ 1. 

To estimate the cumulated deformation in the inner core, we 
assume that the strain rate e found above is applied only on the 
uppermost shear layer of thickness 8. The simplest is to assume 
that both e and 8 are evolving slowly with time, on a timescale 
long compared to the time 8 /M; c over which a layer of thickness $ 
is crystallized. Assuming that the strain rate e(t) is given by u/S 
within the shear layer of thickness 8 and is negligible elsewhere, 
the cumulated deformation is then 


e8 u 
Uic Pe ’ 


(57) 


with u the dimensionless maximum horizontal velocity given by 
Eqs. (48) or (52) depending on the value of the Peclet number. With 
u dimensional, the cumulated strain can also be written as 

(58) 

Uic 


The deformation magnitude below the upper shear layer is given 
by the ratio between the horizontal velocity induced by the Lorentz 
force and the growth rate of the inner core. 

A more elaborated method of estimating e is discussed in 
appendix E. The results are close to what Eq. (57) predicts for 
r > 0.3 ri c (Ti c ), but are more accurate for smaller r. The va¬ 
lidity of both estimates is restricted to conditions under which 
the strong stratification scaling laws applies, which requires that 
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M -C —Ra. This is only verified if the inner core radius is larger 
than ri c (ri c )(Mo/Rao) 1/l ‘ l — 0.02 ri c (ri c ) in the thermal case, 
and ri c (Tic)(Mo/Rcio) 1 / 5 ~ 0.07r; c (T; c ) in the compositional 
case. 


6 APPLICATION TO THE INNER CORE 

To determine in which regime is the inner core, the first step is to 
estimate the ratio —Ra/M 

-Rn = -Apwo = —— ^-3 ( 7.5 X 10 5 . (59) 
M 1 kg.m ' \ Bq J 

Notice that this does not depend on the viscosity. Plausible dimen¬ 
sionless numbers for the Earth’s inner core are obtained from typ¬ 
ical values given in Table 1 and summarized in Table 2. The den¬ 
sity stratification is A p ~ lkg.m -3 irrespectively of the nature of 
the stratification (Deguen & Cardin 2011; Labrosse 2014). Varying 
the parameters within their uncertainty range can change the ratio 
—Ra/M by an order of magnitude at most. The ratio —Ra/M is 
thus unlikely to be smaller than 1, irrespectively of the thermal or 
compositional origin of the stratification. The inner core is strongly 
stratified compared to magnetic forcing. 

If the stratification is of thermal origin, the Peclet number is 
on the order of 1 (Peo = 2.8). Fig. 5 shows that the low Peclet 
number scaling laws still agree reasonably well with the numerical 
results for Peclet numbers around 1; the low Peclet number scal¬ 
ing laws can therefore be used to predict the flow geometry and 
strength in the thermal stratification case. If the stratification is of 
compositional origin, the Peclet number is large (Peo ~ 10 5 ) and 
thus the large Peclet limit scaling laws apply. 

Estimates of the thickness of the upper layer, of the maximum 
velocity in this layer and of the expected strain rate are given in 
Table 3, using values of parameters of Tables 1 and 2. Because the 
viscosity is poorly known, we express these estimates as functions 
of the viscosity. As an exemple, assuming a viscosity of 10 16 Pa.s 
gives a shear layer thickness of 94 km and 60 km in case of ther¬ 
mal and compositional stratification. The velocity in this layer is 
expected to be several orders of magnitude lower than the growth 
rate, and instantaneous deformation due to this flow is small: the 
typical timescale for the deformation is of order 10 2 Gyrs for both 
cases, for rj = 10 16 Pa.s. These values are obtained for the present 
inner core, which means it is the deformation time scale for the 
present uppermost layer. The deformation e ~ u/Pe is a decreas¬ 
ing function of the radius, and thus is higher in depth: compared 
to below the ICB. the strain at r = 0.5 n c is multiplied by 2 for 
thermal stratification, and by 4.6 for compositional stratification. 

As can be seen in Eqs. (55) and (56), the strain rate in the 
shear layer is a decreasing function of the stratification strength. 
This is the opposite of what Deguen et a (20 11 ) found in the case 
of a flow forced by heterogeneous inner core growth (Yoshida et al. 
1996). The flow geometry is similar to what has been found here if 
the inner core is stably stratified, with a shear layer below the ICB 
in which deformation is localized, but, contrary to the case of the 
Lorentz force, the strength of the flow and strain rate increase with 
the strength of the stratification. This difference is due to the fact 
that the velocity is imposed by the boundary conditions at the ICB 
in the case of heterogeneous inner core growth-ease, and therefore 
does not decrease when the stratification strength is increased. In 
contrast, the velocity in the shear layer produced by the Lorentz 
force depends on a balance between the Lorentz force and the vis¬ 
cous forces, and decreases with increasing stratification strength. 


Using the scaling laws developed above (Eqs. (55), (56) and 
(57)), the cumulated strain below the shear layer is given by 


e T ~ 5.6 x 10“ 
and 

ec ~ 2.4 x 10“ 


10 16 Pa.s\ 1/2 


5 , 


3 mT 


10 16 Pa.s \ 2/5 ( B { 


30 


3 mT 


(60) 


(61) 


This shows that a viscosity lower than 10 10 Pa.s is required to ob¬ 
tain a deformation larger than about 1. Such a low viscosity seems 
unrealistic and this suggests that no detectable anisotropy would be 
produced in the bulk of the inner core. 

Using the method in appendix E to estimate the strain, the 
deformation is found to be two orders of magnitude larger close 
to the center of the inner core than at the edge. This means a non- 
negligible strain for viscosity lower than 10 12 Pa.s. 

The uppermost layer has a different behavior because it has- 
does not have enough time to deform. This could stand for an 
isotropic layer at the top of the inner core, as observed by seis¬ 
mic studies. We expect this layer to be of the order of one hundred 
kilometers thick for thermal or compositional stratification. 


7 CONCLUSION AND DISCUSSION 

Following previous studies (Karato 1999; Buffett & Wenk 2001), 
we have developed a complete model for evaluating the deforma¬ 
tion induced by the Lorentz force in a stratified inner core, inves¬ 
tigating the effect of boundary conditions and neutral and strong 
stratification in the case of thermal or compositional stratification. 

Calculating the flow for neutral stratification with different 
mechanical boundary conditions, we show that the boundary con¬ 
ditions depend on the values of the viscosity. If the viscosity is 
low, the ICB acts as an impermeable boundary, with no radial flow 
across the ICB, whereas if the viscosity is large the ICB acts as 
a permeable boundary, with fast melting and solidification at the 
ICB. We find that the velocity is larger than the inner core growth 
rate if the viscosity is lower than 10 lb Pa.s. Unlike previous studies, 
the boundary conditions assumed here are of impermeable type. 

If the inner core has a stable density stratification, then we 
find that the stratification strongly alters the flow induced by the 
poloidal component of the Lorentz force. The deformation is con¬ 
centrated in a thin shear layer at the top of the inner core, which 
thickness does not depend on the magnetic field strength, but de¬ 
pends on both the density stratification and the Peclet number, 
which compares the timescales of inner core growth and diffusion. 

However, the deformation rate in this regime is predicted to be 
too small for producing significant LPO in most of the inner core, 
unless the inner core viscosity is smaller than lO 10 — 10 12 Pa.s. 
The cumulated deformation can be two orders of magnitude larger 
close to the center of the inner core, but remains smaller than 1 if 
the inner core viscosity is larger than 10 12 Pa.s. 

We have made a number of simplifying assumption, but relax¬ 
ing them is unlikely to significantly alter our conclusions. Our es¬ 
timated values of the deformation are probably upper-bounds. The 
effective strain in the inner core induced by the poloidal Lorentz 
force is expected to be even smaller than these values. Indeed, we 
use assumptions that maximize the strain rate in the inner core, (i) 
The geometry and strength of the magnetic field have been chosen 
to maximize the effect of the Lorentz force: the degree (2,0) pen¬ 
etrates deeper in the inner core than higher orders (smaller length- 
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Table 3. Estimates of the thickness, maximal horizontal velocity and strain rate of the upper layer for thermal stratification (low Pe) and compositional 
stratification (large Pe). 


Thermal stratification, low Peclet 


Compositional stratification, large Peclet 


Thickness S 

Maximal horizontal velocity 0 u 
Instantaneous strain rate b e 
Strain e = u/Pe 


rj \ 1/6 


' 10 16 Pa.s 


X 94 km 


' 10 16 Pa.s\ 1/2 


X 2.2 X 10 —14 m.s —1 


f 10 16 Pa.s\ 2/3 , 

I -) x 7.4 x 10” 12 yrs^ 1 

/ 10 16 Pa.s \ 1/2 r ic , 

- — X 5.6 X 1CT 4 


ri \ V s 

—-z -) X 60 km 

10 16 Pa.s/ 

10 16 Pa.s\ 2/5 ,, , 

- x 0.9 x 10~ 14 m.s -1 

J 

10 16 Pa.s\ 3/5 


10 1 ® Pa.s 
V 


x 4.8 x 10~ 12 yrs - 


2/5 


11/5 


X 2.4 X 10’ 


this value has to be compared with a typical value for the growth rate: Ui c (ri c ) Ri 10 11 m.s 


’at t = tv. 


scales) and is less likely to vary with time. The smaller scale com¬ 
ponents of the magnetic field vary on a shorter timescale, which 
make them more sensitive to skin effect, (ii) We have assumed the 
magnetic field to be time-independent, a reasonable assumption be¬ 
cause the fluctuations associated with outer core dynamics occur on 
a timescale short compared to the timescale of inner core dynam¬ 
ics. The evolution of the magnetic field strength on the timescale 
of inner core growth is poorly constrained, but seems unlikely to 
have involved order of magnitude variations, (iii) Though a growth 
law of the form ri c oc t 1 ' 2 has been assumed in the numerical 
simulations, our derivation of the scaling laws (section 4.2) makes 
no assumption on the inner core growth law, and our scaling law 
should therefore also apply if this assumption is relaxed. 

Finally, we have focused on the effect of the azimuthal compo¬ 
nent of the magnetic field, which produces a poloidal Lorentz force 
in the inner core, and have left aside the combined effect of the 
azimuthal and 2 —component (parallel to Earth's spin axis), which 
produce an azimuthal Lorentz force driving an azimuthal flow (Buf¬ 
fett & Wenk 2001 ). There is no loss of generality involved, because 
the axisymmetric poloidal flow we have investigated and the az¬ 
imuthal flow calculated by Buffett & Wenk (2001) are perfectly de¬ 
coupled, and add up linearly. The azimuthal flow is horizontal, and 
is therefore not affected by the thermal and compositional fields and 
their perturbations by the axisymmetric poloidal flow. Conversely, 
since the flow and density perturbations induced by the axisymmet¬ 
ric poloidal Lorentz force are axisymmetric, they are not affected 
by an azimuthal flow. The azimuthal flow velocity is 


1 B Z B$ r 3 

10 Mo V tfc 


sin 0 


(62) 


(Buffett & Wenk 2001), and the associated strain rate is 


e(r, 9) = 9) = BzB ^ s i n Q. (63) 

10 /iotMTc 

This is always larger than the strain rate predicted by our scaling 
laws in the strongly stratified regime (Eqs. (55) and (56) for ther¬ 
mal and compositional stratification), which implies that this de¬ 
formation field will dominate over deformation due to the poloidal 
component of the Lorentz force. 
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APPENDIX A: POLOIDAL/TOROIDAL 
DECOMPOSITION 

A1 Momentum equation using poloidal decomposition 

Following Ricard & Vigny (1989) and Ribe (2007), a flow induced 
by internal density anomalies in a constant viscosity fluid is purely 
poloidal in a spherical shell, if the surface boundary conditions also 
have a zero vertical vorticity. Considering also that the external 
forcing by the Lorentz force is purely poloidal (r ■ V x Fl = 0), 
no toroidal flow is expected in our problem. The velocity field can 
thus be written by introducing the scalar P defined such that 

u = V x V x ( Pr ), (A.l) 


where r = re r is the position vector. 

Applying r • (V x Vx ) to Eq. (13) , we obtain that 

0 = Ra(t)L 2 Q - (V 2 ) 2 L 2 P + M(t)r ■ (V x V x F L ), (A.2) 

where L 2 is the operator defined by 

^ a =--44 ( sin ^4) (A3) 

sin 9 o9 \ 39 ) sin 9 d<p 2 

The last term of Eq. (A.2) is computed from the Lorentz force 
defined in Eq. (3) and gives 

r ■ (V x V x F l ) = 8r 2 (l - 3 cos 2 9) = ~^= r 2 Yn°, (A.4) 

v 5 


where Y§ = ^ (3 cos 2 61-1). 
Eq. (A.2) becomes 


0 = Ra(t)L 2 Q - (V 2 ) 2 L 2 P- M(t)^r 2 Y° 


(A.5) 


When expanding the two scalar field 0 and P with spherical 
harmonics yj m (0, (j>) that satisfy L 2 V; m = 1(1 +1 )Y{ n , we defined 
new variables if 1 and p™ by 


0 = tr(r)Yr, (A. 6 ) 

P = pT{r)Yi m , l> 1. (A.7) 


Eq. (A.5) is eventually written as 


In term of poloidal decomposition of the velocity field, the set 
of equations for the boundary conditions at r = 1 is modified as 


rl 2 n m 

+ [!(! + !)- 2 ^= 0 , l> 1 , 


(A. 13) 


dr 3 


31(1 + 1) 


1 dpR 

r dr 


l(l + l)V(t) 


6 


m 


Pi , 


l > l. 

(A. 14) 


APPENDIX B: THERMAE STRATIFICATION 

We derive here the expression of the reference diffusive poten¬ 
tial temperature profile given in Eq. (24), under the assumption of 

Tic oc t 1/2 . 

The reference potential temperature is the solution of 

dO , u ic (t) dQ _ k t Id / 2 dO\ q , R . 

dt V r ic (t) dr rf 0 (f) dr \ dr) T ’ 

obtained from Eq. (11) by taking \ = 0 an d u = 0. The source 
term St is constant if ri c oc t 1 ' 2 . The potential temperature 0 
is a function of r, kt , ri c (t), St and m c (t ) only. According to 
the Vaschy-Buckingham theorem, we can form only three inde¬ 
pendent dimensionless groups from these variables (6 dimensional 
variables - 3 independent physical units), one possible set being 
0 / (STrfj kt), UicUc/KT . and r. With rt c oc t 1 ' 2 , Ui c ri c /KT is 
constant and equal to Pe o. The potential temperature must there¬ 
fore be of the form 

Q = S T r 2 c (i) /(f ,p eTo) (B.2) 

By definition, the potential temperature is equal to 0 at the 
ICB, which implies f(f = l,Pero) = 0. With ri c (t) = 
ric(Tic)(t/T ic ) 1/ 2 and noting that r 2 c (r ic )/ (K T T ic ) = £to = 
2PeTQ (see Eqs. (17) and (21)), inserting Eq. (B.2) into Eq. (B.l) 
yields 

0 = f" + /' 0 + Pe T0 ^j - 2Pe TO f + 1, (B.3) 


DhT + ^ M(t)r 2 8 2l 5 0m 

V 5/(t + 1) 


Ra(t)t™ = 0, l > 1, 

(A. 8 ) 


where S is the Kronecker symbol and Di is the second order differ¬ 
ential operator defined by 


d 2 2d 
dr 2 r dr 


1(1 + 1 ) 


(A.9) 


where f' and f" stand for the first and second derivatives of / 
with respect to r. Looking for a polynomial solution in r satisfying 
f(r = 1, Pe to) = 0 , we find that / = (1 — f 2 )/( 6 + 2Pero), 
which gives 

= S T r 2 c (t) 

6kt(1 + Pero/3) 



A2 Poloidal decomposition of the boundary conditions 


From Deguen et al. (2013), the boundary conditions at r = 1 are 
written as 




d 

(ue 

\ , l 

du r 

= 0, 

TrB 

= V 

r dr 


) + r 

~d9 


r d 

(u$\ 


1 

du r 

= 0, 

= V 

r dr 

(v) 

+ 

rsind 

d(j> 

— V(t)(u r — U 

ic) 

du, 

dr 

+ p 

' = 0 . 


(A. 10) 
(A.l 1) 
(A. 12) 


P(t) is the dimensionless parameter that characterizes the resis¬ 
tance to phase change as defined in Eq. (35). 


APPENDIX C: COMPOSITIONAL STRATIFICATION 

The source term S c of the conservation of light elements is directly 
related to the evolution of the concentration of light elements in the 
solid that freezes at the inner core boundary cf c (f). Following Gub- 
bins et al. (2013) and Labrosse (2014), this term depends both on 
the evolution of the concentration in the liquid outer core, which 
increases when the inner core grows because the solid incorporates 
less light elements than is present in the outer core, and on the evo¬ 
lution of the partition coefficient between solid and liquid. 

In this paper, we will focus on the simplest case for which the 
partition coefficient does not depend on temperature or concentra¬ 
tion. Thus, the concentration in the solid is increasing with the ra¬ 
dius of the inner core, as the concentration in the liquid increases. 
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This will promote a stably stratified inner core, whereas Gubbins 
et al. (2013) and Labrosse (2014) focused on the potentially desta¬ 
bilizing effects of variations of the partition coefficient. 

To estimate the light element concentration, we note M c = 
Mi c + M oc the total mass of the Earth’s core. When increasing 
the inner core mass by d Mi c , the mass of the outer core light el¬ 
ements decreases by d(c l M oc ) = —cf c dMi c . The total mass of 
the Earth’s core is constant, which gives d M oc = — d Mi c and 


d c l 
c l 


(k - 1 ) 


dM 0 , 

M 0 p. 


(C.l) 


where k is the partition coefficient defined as k = cfo/c 1 . 

Eq. (C.l) can be integrated with the assumption of a constant 
partition coefficient. Integration between ( cb,M c ) and (c , M oc ), 
corresponding to before the inner core formation and any time after, 
this gives 


c l (t) = Co 



(C.2) 


When ignoring radial density variations in the outer core, the 
ratio M oc /M c is simply 1 — ( ri c /r c ) 3 . Taking into account com¬ 
pressibility (radial density variations in the core) and the density 
jump at the ICB results in a stratification approximately 15% larger. 
The light element concentration at the inner core boundary is thus 
directly obtained from the liquid concentration as 


C; c (i) = kc 0 



fjc(t) 

r c 



(C.3) 


Table Al. Extreme values for RMS velocity and maximum of the absolute 
value of the horizontal velocity for two extreme values of V. Velocities are 
proportional to M and thus only the value v/M is given. 



V 0 

V —>• 00 

Vrms / M 

0.06609 

0.00805 

max | ug(r, 0)\/M 

0.06944 

0.01270 


From the coefficients p ™, the vertical and horizontal velocities 
are 

m 

Ur = Y,l(l + l)^Yr, (D. 6 ) 

l,m 



(D.7) 


with YC the surface spherical harmonics. 

The root mean square velocity ( V rrns ) of the system is defined 
as 

o /* 2 (/> pit /*1 

Vrms — ~r~ / / / {if.+ u 2 g) sin Or 2 dr ddd(j>. (D. 8 ) 

4 * Jo Jo Jo 

From Eq. (D.7) and (D. 8 ), we obtain the maximum absolute 
value of the horizontal velocity and the RMS velocity that are 
shown on Fig. 2. Both graphs have a sigmoid shape, and thus we 
are mostly interested in the extreme values for each velocity, which 
are given in Table Al. 


APPENDIX D: ANALYTIC SOLUTION FOR NEUTRAL 
STRATIFICATION 


We solve Eq. (A. 8 ) for a neutral stratification, Ra = 0, with the 
boundary conditions (A. 13) and (A. 14) described in section (A2). 
Equation (A. 8 ) is thus written as 

Dtp? + ™ Mr* 8216 on = 0 , (D.l) 

V5 1 (l + 1) 

and can be solved analytically. 

Except for (l = 2 ,m = 0), p™ = 0 is solution of the 
Eq. (D.l) and verifies the boundary conditions (A.13) and (A.14). 
For (l = 2, m = 0), we have 




(D.2) 


and for r = 1 


d 2 p ° 2 
dr 2 


+4 


r «_ 18 l iA 

dr 3 r dr 



(D.3) 

(D.4) 


Eq. (D.2) is solved considering a sum of polynomial functions, 
and adding the boundary conditions (D.3) and (D.4), we obtain the 
coefficient p% as 


APPENDIX E: INTEGRATION OVER TIME OF THE 
DEFORMATION 

El General discussion 

In general, the texturation mechanism is a nonlinear process, but an 
upper bound of the total deformation can be inferred by considering 
that the strain adds up linearly. The material is deformed at a strain 
rate e during the time 5/ui c needed to grow a layer of thickness 5, 
and so 

e(rW)= ^ ) S) = A’ (R1) 

with Pe = Ui c ric/ ft. 

This equation leads to simple forms at the low and large Peclet 
limits, with e oc t -1 ' 2 for thermal stratification and low Peclet, and 
e oc i -11 / 10 for compositional stratification and large Peclet. In 
what follows, we will compare the simple estimate given above 
with results of more elaborate calculations. 

The total deformation of a given material during a time r 
can be inferred more precisely by J () r k(t)dt. Using dimensionless 
quantities described in the main sections, the deformation of a strat¬ 
ified sphere subject to a magnetic forcing is 

e(r) = [ e(r, t ) dt, (E.2) 

Jo 


P°(r) = 


M 


3 3 745 

14 4 

+ ~T r 
5 


9 ^ 2 ( 204 r 4 
5 r + IT 19 + 5P 


544 r 2 \ 
19 + 5 PJ ' 


(D.5) 


with e the strain rate function that will be described by a rectangular 
function as 


e(r, f) = 


e vM{t) r 2 

0 


if ri c (f)(l — 8 ) < r < r ic (t ) 

elsewhere. 
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(E.3) 


The estimations of e V M depend on the scaling laws defined in 
section 5, and also on the time dependence of the parameters we 
have defined. 

Because 7y c (t)( 1 — 8 ) < ri c (t) Vf, integrating over time 
the function defined by (E.3) is equivalent to integrate it be¬ 
tween tmin(r) and t m ax(r), where t mi „ and tmax are defined by 
ric(tmax){ 1 — 8) = r and — r, 



£vM (f ) 


rl(t) 


dt. 


(E.4) 


E2 Low Pe - Thermal stratification 

In the low Peclet limit, the dimensionless thickness, maximal hori¬ 
zontal velocity and strain rate of the uppermost layer are given by 

8 ~ {-Ra)~ 1/6 , (E.5) 

M(-Ra)~ 1/2 , (E.6) 

e ~ M(-Ra)~ 1/3 . (E.7) 

In dimensional form, 5 happens to be constant with time for 
thermal stratification 


8 = 1.9643ri c (Ti c )(— Rap) 1/6 . 


Thus, tmi n and t max are easy to defined as 

• {r) = T “(^kj' 


(E.8) 


(E.9) 


tmax (r) = T ic ( T ^ ^ , (E.10) 

\ric{Tic)) 

except for time close to Ti c because the inner core has not enough 
time to be deformed, and t max = Ti c . For small radius, the limit 
will be defined by 8 = 1, which is here —Ra(t) = M(t), about 
28km for typical values of the parameters. 

For thermal stratification and time dependence as defined pre¬ 
viously, the instantaneous deformation is 


e(t)—^m = 0.2148M 0 (—flo 0 )~ 1/3 * n of" 1 , (E.ll) 

r ic[t) r ic\ T ic) 

= £o Tic t 1 , (E. 12) 

with £o in s _1 and t in s. 


f e 0 T ic 2 1 n^, 

e(r) 'i 


for r < ri c (ri c ) - 5, 
for r > ndjic) — 8, 


(E.13) 


The strain rate is assumed to be constant over the layer 8 
whereas we could have used numerical results of the simulations 
to have the exact repartition and profile of strain over radius and 
time. But because the von Mises strain rate profile over radius is 
close to a rectangular function, it is easier to work with an analytic 
solution such as the one discussed here. It implies a linear increase 
of the absolute value of the strain, which is unlikely. This will in 
general overestimate the total strain. 

Comparison between the strain computed from (E.13) and the 
simplest solution u/Pe discussed in the text is plotted on Fig. Al. 
The solution u/Pe is a good approximation for radius larger than 
0.3 except in the uppermost layer, in which the deformation did 
not occur completely yet. It is interesting to notice that this mag¬ 
netic forcing is expected to be several orders of magnitude more 
efficient when the inner core was younger. 



Figure Al. Strain as a function of the radius of the inner core for ther¬ 
mal stratification (a.) and compositional stratification (b.). Integration from 
equation (E.4) is the red line (analytical solution for thermal stratification 
and numerical solution for compositional stratification), and the blue lines 
stand for the estimation u/Pe, which is valid for r > 0.5n c . The min¬ 
imum radius is computed for Ra/M = 1, limit under which the strong 
stratification approximation is no longer valid. 


E3 Large Pe - Compositional stratification 


For compositional stratification, 



1/5 

ric(t"ic) 



(E. 14) 


No exact solution for inverting r ic (r ic )(t m ax/Tic) 1 ^ 2 — 
<5(Ti c )(t/Tic) -1 / 5 = r can be found. Fig. Al b. shows the strain 
rate according to numerical integration and the approximation 
u/Pe. u/Pe is a good approximation for r > 0.3 n c . 















